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PACS 03 . 75 . Kk - Dynamic properties of condensates; collective and hydrodynamic excitations, 
superfluid flow 

PACS 71.10.Fd - Lattice fermion models (Hubbard model, etc.) 

Abstract. - We use the density-matrix renormalization group method to investigate ground-state 
and dynamic properties of the one-dimensional Bose-Hubbard model, the effective model of ultra- 
cold bosonic atoms in an optical lattice. For fixed maximum site occupancy rij, = 5, we calculate 
the phase boundaries between the Mott insulator and the 'superfluid' phase for the lowest two 
Mott lobes. We extract the Tomonaga-Luttinger parameter from the density-density correlation 
function and determine accurately the critical interaction strength for the Mott transition. For 
both phases, we study the momentum distribution function in the homogeneous system, and the 
particle distribution and quasi-momentum distribution functions in a parabolic trap. With our 
zero-temperature method we determine the photoemission spectra in the Mott insulator and in 
the 'superfluid' phase of the one-dimensional Bose-Hubbard model. In the insulator, the Mott 
gap separates the quasi-particle and quasi-hole dispersions. In the 'superfluid' phase the spectral 
weight is concentrated around zero momentum. 
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Introduction. — At very low temperatures, bosonic 
atoms which are loaded into an optical lattice become su- 
perfluid for a shallow optical potential and Mott insula- 
tors for a deep optical potential. The transition between 
both phases has been observed experimentally [TJ ; for a re- 
cent review, see [2J. The Bose-Hubbard model provides a 
reasonable description of the experimental situation, and 
its ground-state phase diagram in two and three dimen- 
sions has been determined fairly accurately by perturba- 
tion theory [3]-[5] and quantum Monte-Carlo (QMC) cal- 
culations [6]-[8]. 

Bosons on a chain are also accessible experimentally [0. 
so that it is interesting to study the one-dimensional Bose- 
Hubbard model. The physics in one dimension is rather 
peculiar. For example, the state with the lowest kinetic en- 
ergy is not macroscopically occupied in the 'superfluid' 111 
but it is characterised by an algebraic divergence of the 
momentum distribution; for a review, see [11] , Moreover, 
the Mott gap is exponentially small in the Mott insula- 
tor close to the phase transition. Therefore, it is very 
difficult to determine the critical interaction strength nu- 
merically. This problem also impairs the applicability of 
strong-coupling perturbation theory. 

fn one dimension and at zero temperature, the density- 
matrix renormalisation group (DMRG) method [T2HT¥| 



permits the calculation of ground-state properties with an 
excellent accuracy for large systems so that the extrapo- 
lation to the thermodynamic limit can be performed reli- 
ably. In this work, we use the density-density correlation 
function to calculate the Tomonaga-Luttinger parameter 
from which we determine the Mott transition accurately. 
Moreover, we obtain the momentum distribution and the 
particle distribution for bosons on a homogeneous chain 
and in the presence of a harmonic trap. Using the dynam- 
ical DMRG [15], we calculate the single-particle spectral 
function at zero temperature in the 'superfluid' and the 
Mott insulating phases. 

Bose-Hubbard model. — The Hamilton operator 
for the Bose-Hubbard model on a chain with an even num- 
ber of sites L in a harmonic potential is defined by 

# = -tEf^+i+^+AO + lE^- 1 ) 



+kEc?-o 2 



(i) 



where St and 6 • are the creation and annihilation operators 

for bosons on site j, fij = fetfc is the boson number opera- 
tor on site j, t is the tunnel amplitude between neighbour- 
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ing lattice sites, U > denotes the strength of the on-site 
Coulomb repulsion, V c parameterises the curvature of the 
quadratic confining potential, and r c = (L + l)/2 denotes 
the central position of the chain. In the following, we set 
U = 1 as our energy unit, unless stated otherwise. 

Constrained Bose-Hubbard model. In general, the 
Bosc- Hubbard model cannot be solved analytically. In the 
low-density limit, the model reduces to the bose gas with 
(5-potential interaction which was solved by Lieb and Lin- 
iger |16| . The fact that three or more bosons may occupy 
the same site forms the major obstacle on the way to an ex- 
act solution. Since multiple-occupancies also pose techni- 
cal problems in numerical approaches, the Bose-Hubbard 
model is usually approximated by the constraint that there 
is a maximal number of bosons per site, < nb < J\f — 1. 
This constrained Bose-Hubbard model has TV degrees of 
freedom per site so that it can be written in terms of spin 
variables with S = (TV - l)/2. The case TV = 2 is triv- 
ial because the hard-core Bose-Hubbard model has no in- 
teraction term. It reduces to a model for free spinless 
fermions whose properties are known exactly j!7) . The 
Bose-Hubbard model is recovered in the limit TV — > oo. 
In general, however, the SU(TV)-Bethe Ansatz equations 
do not solve the constrained Bose-Hubbard model [TSlfT^] . 

In our work, we study the restricted Bose-Hubbard 
model with TV = 6, i.e., rib < 5. Our results are repre- 
sentative for the original Bose-Hubbard model ([T]) because 
multiple lattice occupancies are strongly suppressed in the 
parameter regions of interest to us, U/t > 2 and fillings 
p = N/L < n b . 

Numerical algorithm. We adopt the DMRG 
method [TJ] as our numerical tool for the calcula- 
tion of ground-state properties for constrained bose 
systems [T3J[l4]. For the spectral properties, we employ 
the dynamical DMRG (DDMRG) [15]. 

The considered lattices are large enough to permit reli- 
able extrapolations to the thermodynamic limit for the 
physical quantities of interest to us. We keep up to 
m = 2000 density-matrix eigenstates, so that the dis- 
carded weight is always smaller than 1 x 10 -10 . 

We checked our algorithm for rib = 1 against the exact 
result [T7] . The exact ground-state energy in the thermo- 
dynamical limit and the extrapolated ground-state energy 
from DMRG agree to four-digit accuracy. 

Ground-state phase diagram. At integer filling p = 
N/L, the Bose-Hubbard model in one dimension describes 
a Mott transition between the 'superfluid' phase, char- 
acterised by a divergence of the momentum distribution 
at momentum k = [TP] , and a Mott insulating phase, 
characterised by a finite gap for single-particle excitations. 
The latter is defined by the energy difference between the 
chemical potentials for half band filling and one particle 
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Fig. 1: Phase diagram of the one-dimensional constrained 
Bose-Hubbard model (nj < 5) from DMRG with 'superfluid' 
(SF) and Mott insulating (MI) regions. The symbols confine 
the regions with a finite Mott gap, A > 0, extrapolated from 
A(L) for L < 128. The position of the Mott tips has been 
obtained from the Tomonaga-Luttinger parameter. 

less than half filling, 

A(L) = p + (L)-p-(L), 
p+(L) = E (L,N+1)-E (L,N) , (2) 
pT{L) = E (L,N)-E (L,N-1) , 

where Eq(L, TV) is the ground-state energy for L sites and 
TV particles. In the thermodynamical limit, TV, L — > oo and 
p = N/L integer, the gap is finite for the Mott insulator, 
A = lim^r .l^oo > 0, so that the system becomes 

incompressible when we go from the 'superfluid' phase to 
the Mott insulating phase. 

The Mott transition lines in the p-U ground-state 
phase diagram have been previously determined by various 
analytical and numerical methods, e.g., strong-coupling 
expansions [30] HI] > variational cluster approach [22] , 
QMC [23l[24], and DMRG [13j[T4j. In fig. Q] we show the 
phase diagram for the first Mott lobe (p = 1) and the 
second Mott lobe (p = 2) as obtained from our DMRG 
calculations with system sizes up to L = 128. 

The overall shape of the Mott lobes agrees with pre- 
vious results. Here, we provide accurate data for the 
second Mott lobe, and the values for the critical interac- 
tion strength for the first two Mott lobes which we obtain 
from the Tomonaga-Luttinger parameter. At the tip of 
each Mott lobe, the model is in the universality class of 
the XY spin model so that there is a Kosterlitz-Thouless 
phase transition with the Tomonaga-Luttinger parameter 
Kb =1/2, and the gap is exponentially small in the vicin- 
ity of (t/U) c . In contrast, SU(TV)-Bethe Ansatz equations 
predict a discontinuity of the gap at the critical interaction 

for tv > 3 nam]. 
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Ground-state properties. 

Tomonaga-Luttinger parameter and critical interac- 
tions for the Mott transition. The low-energy excita- 
tions of interacting bosons in the superfluid phase are 
gapless linear excitations ('phonons'). As in the case of 
fcrmionic systems in one dimension J25[[26], the Tomona- 
ga-Luttinger parameter Kb determines the asymptotic be- 
haviour of the correlation functions in the 'superfluid' 
phase, and various correlations functions have been used 
to extract Kb [2"Tl430j . Here, we employ the density-density 
correlation function which is defined by the ground-state 
expectation value 



1 L 

C(r) = —y](ht +r fit) - (n e+r )(n e ) . 
1=1 

Asymptotically, it behaves like 

1 1 Ap 2 cos(2npr) 



C(r ->■ oo) 



2K b (7rr) : 



(3) 



(4) 



Thus, we can extract Kb from the derivative of its Fourier 
transformation, 

L 

C(q)=Y / e- iqr C(r),0<q<2n ) (5) 

r=l 

as q = 0. In the thermodynamic limit one finds 

-J- = l im ^M. (6) 

2ttRi, <?->o q 

In order to treat finite systems in numerical calcula- 
tions [25] . we translate (JSJ into 



, r . = lim —G 
2K b (L) l^oo 2 



(7) 



and extrapolate K\,(V) to the thermodynamical limit. 

In refs. [H] and [3T] the transition point has been also 
determined from the Luttinger parameter Kb- However, 
these authors estimated Kb from the single-particle den- 
sity matrix 



r(r) = $%) ~ r ~ Kb/2 for r > 1 



(8) 



In their work, the extrapolation for the critical point t c 
depends on the interval used for the fits to T(r), see table I 
in ref. [14j . When we derive the Luttinger parameter from 
cq. © we can avoid this problem. 

As shown in fig. [21 Kb{L) can be reliably extrapolated 
to the thermodynamic limit using polynomial functions in 
1/L. For p = 1, we clearly have K b (t/U = 0.3) > 1/2. 
When we extrapolate our data for up to L = 1024 lat- 
tice sites, we find K b (t/U = 0.304) > 1/2 but K b (t/U = 
0.306) < 1/2. Therefore, we locate the transition point at 
t c = 0.305 ±0.001 for the first Mott lobe. In the same way 
we find the transition point for the second Mott lobe at 
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Fig. 2: Finite-size scaling for the Tomonaga-Luttinger pa- 
rameter Kb in the one-dimensional constrained Bose-Hubbard 
model (rib = 5) for the first (p = 1) and second (p = 2) Mott 
lobes, using the DMRG with open boundary conditions. The 
lines are polynomial fits. The insets give the n(,-dependence of 
C(2ir/L) for various system sizes at t/U — 0.305 (left panel) 
and t/U = 0.18 (right panel). 



t c = 0.180 ±0.001 for the restricted Bose-Hubbard model 
with nb < 5. Using the same method, we have verified nu- 
merically that the values for the critical coupling (t/U) c 
are the same for nb = 4, 6 within our extrapolation uncer- 
tainty. 

Note that Kb(t < t c ) is not defined because we are in the 
Mott insulating phase. However, Kb(L) is finite and con- 
tinuous over the Kosterlitz-Thouless transition because 
the Mott gap is exponentially small near t c . Nevertheless, 
our approach remains applicable as has been shown for 
various fermionic models in refs. [2"rJll32j . 

Previous groups located the Kosterlitz-Thouless transi- 
tion for the first Mott lobe at values consistent with ours. 
In their DMRG work jl4] . Kiihner et al. computed the 
Luttinger parameter using their DMRG algorithm on lat- 
tices with up to L — 1024 sites. From their fit to T(r), 
cq. ©, they found t c = 0.297 ± 0.01. Based on the same 
correlation function, Zakrzewski and Delande [3T] gave 
t c = 0.2975 ± 0.005 for the first and t c = 0.175 ± 0.002 
for the second Mott lobe for nb — 6. The determined de- 
values of such a kind significantly depend on the interval 
of r, which is not the case within our approach. Lauchli 
and Kollath J3HJ determined the critical point from the 
block entropy, by combining the recently developed quan- 
tum information theory with the DMRG. Our result is 
within their region of the estimated values for t c (see fig. 2 
in ref. [33] L In a combination of an exact diagonalisation 
study for systems with up to L — 12 sites and a renormali- 
sation group approach, Kashurnikov and Svistunov found 
t c = 0.304 ± 0.002 [34], and their QMC calculations to- 
gether with Kravasin gave t c = 0.300±0.005 [23]- Another 
QMC calculation in combination with a renormalisation- 
group flow analysis of the finite-temperature data gave 
t c = 0.305(4) [35], in perfect agreement with the result of 
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Fig. 3: (Colour online) Finite-size dependence of the momen- 
tum distribution function n(k) for t = 0.1 (Mott insulator) and 
t = 0.305 ('superfluid') in the one-dimensional Bose-Hubbard 
model using the DMRG with periodic boundary conditions. 
Inset: n(k) for L — 64. The solid line in the inset gives the 
result of strong-coupling theory to third order, eq. (|10|) [5ll37|. 



Fig. 4: (Colour online) Occupation probabilities in the one- 
dimensional constrained Bose-Hubbard model (rib < 5) in a 
parabolic trap potential of strength V c /U = 0.001. We show 
the results for N = 24, 40, 54 and L = 100 (p = 0.24, 0.40, 0.54) 
for t/U = 0.125 for (a) the local densities and (b) the 

pseudo- momentum distribution n(q). 



our zero-temperature study. 

Momentum distribution function. Using the DMRG, 
the momentum distribution function n(k) can be calcu- 
lated by taking the Fourier transformation of the single- 
particle density matrix 

j,£=l 

where k = 2-Km/L for m = —L/2 — 1, . . . , L/2 holds for pe- 
riodic boundary conditions [36] . Note that the momentum 
distribution function fulfils the sum rule n(k) — N. In 
all cases of fig.[3]the numerical deviation £ = \N— J^k n (k)\ 
is always small, £ < 1.0 x 10~ 3 . 

The difference between the superfluid phase and the 
Mott insulator is most markedly seen in the momentum 
distribution n(k) at momentum k = 0: in the insulating 
phase, n(k = 0) remains finite whereas it diverges as a 
function of system size in the superfluid phase, as shown 
in fig. El At t/U = 0.1, n(k = 0) is almost independent 
of system size, and the momentum distribution n(k) is a 
smooth function of momentum k. Strong-coupling pertur- 
bation theory to third order 0137] predicts (x = t/U) 

n [3] (k) = 1 + 2Ci cos(fc) + 2C* 2 cos(2fc) + 2C* 3 cos(3fc) , 
Ci = Ax - 8x 3 , C 2 = 18a; 2 , C 3 = 88a; 3 . (10) 

Our numerical results for t/U = 0.1 favourably compare 
with this expression, see the inset of fig. [3] 

At t/U = 0.305, above the critical point, n(k = 0) in- 
creases rapidly with system size. In one spatial dimension 
there is no true superfluid with a macroscopic value for 
n(k = 0) in the thermodynamic limit |10j . Instead, we 
have from © n(\k\ -> 0) - \k\~ u , v=\- K b /2 < 1. 

Local densities for the Bose-Hubbard model in a trap. 
In the presence of the confining potential V c in the 



model ([T]) , the density profile over the trap is no longer ho- 
mogeneous, see, e.g., ref. [38]. For an open system, we de- 
fine the quasi-momentum distribution n(q) = (w(q)b(q)) 
using the quasi-momentum states of particles in a box, 

£ 

with q = Trn q /(L + 1) for integers 1 < n q < L. 

As demonstrated by Batrouni et al. |38j and Kollath et 
al. |36j , the potential confines the particles in the middle of 
the trap. For small fillings, the local occupancies display 
a bell-shaped distribution, where the maximum does not 
reach the Mott plateau value, {b)bj){p = 0.24) < 1. The 
quasi-momentum distribution n(q) for this superfluid in a 
trap shows a prominent peak at k = 0. For a larger filling, 
p = 0.40, there exists a Mott plateau, (£>jS,-)(p = 0.40) = 1 
for 40 < j < 60. Recall that, for t/U = 0.125, the homo- 
geneous system at filling p = 1 is a Mott insulator. Cor- 
respondingly, the peak in the pseudo-momentum distribu- 
tion at k = is smaller for p = 0.40 than for p = 0.24, 
see fig. S] (b). Finally, at filling p = 0.54, the confin- 
ing potential and the bosons' tendency to cluster over- 
come the repulsive potential in the middle of the trap so 
that local occupancies larger than unity arc seen inside the 
trap. Correspondingly, the peak intensity of the pseudo- 
momentum distribution at p = 0.54 exceeds its value for 
p = 0.40. 

Photoemission spectra. — Single-particle excita- 
tions associated with the injection or emission of a boson 
with wave vector q and frequency ui, A + (q,u>) or A~(q, ui), 
are described by the spectral functions 

A ± (q,u 1 ) = £ K^lS^g)!^) 2 *^ T u>±) , (12) 

n 

where b + (q) = b'(q) and b (q) = b(q) create/annihilate 
particles with pseudo-momentum q. Moreover, \ipo) is the 
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Fig. 5: Intensity (left panels) and line-shape (right panels) of 
the single-boson spectral functions A(q, uf) in the Mott insulat- 
ing (MI) phase for t/U = 0.05 (upper panels) and t/U = 0.1 
(lower panels) with system size L = 64 at filling p = 1 using 
the DDMRG technique with open boundary conditions. 

ground state of a L-site system in the iV-particle sector 
while |^) denote the nth excited states in the (TV ± 1)- 
particle sectors with excitation energies = — E Q . 

So far, very few data are available for the (inverse) pho- 
toemission spectra in the one-dimensional (constrained) 
Bose-Hubbard model. Analytical results include the vari- 
ational cluster perturbation theory [22] , the random phase 
approximation [39) , and strong-coupling theory [40[ . Pip- 
pan et al. [H] combined QMC at low but finite temper- 
atures with the maximum-entropy method to extract the 
spectral functions. 

In the following we present the (inverse)-photoemission 
spectra at zero temperature using the numerically exact 
dynamical DMRG method [HUH]. We keep m = 500 
states to obtain the ground state in the first five DMRG 
sweeps and take m = 200 states for the calculation of 
the various spectra from ([P2"| by DDMRG. For a bosonic 
system the following sum-rules hold, 

/oo 
duj(A + (k,uj) - A-(k,u)) = 1, (13) 
OO 




dw (A + {k,u) + A-(k,u)) = n(k) . (14) 



In our DDMRG calculations, both sum-rules are fulfilled 
with high precision. 

In fig. [5] we show the results for the Mott insulator with 
p = 1. The spectra A(q, u>) = A + (q, uj) + A~(q, ui) for fixed 
q consist of two Lorentzians of width rj = 0.04, the size 
of the broadening introduced in the DDMRG procedure. 
The quality of the fits suggests that the quasi-particle life- 
time is very large in the Mott insulator. 
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Fig. 6: (Colour online) Quasi-particle dispersions u(k) in the 
Mott insulating phase at filling p — 1 for t/U = 0.05 and 
t/U — 0.1 from Lorentz-fits to the spectral functions A(q,u>). 
For comparison, we also show the strong-coupling dispersions 
for the propagation of a hole and a double occupancy, cjh. P (k), 
and the mean-field result of ref. [ID]. Inset: A(q,u>) for t/U = 
0.1 at q = 7r/33 with L = 32 (line) and q = 2n/65 with L = 64 
(circles) demonstrating the negligible system size dependence. 

In fig.[5]we show the quasi-particle dispersions for t/U = 
0.05 and t/U = 0.1 which we extracted from the fits of the 
spectral functions to two Lorentz peaks at ui = oj^(k). For 
comparison, we include the mean-field result [40] and the 
strong-coupling result. For large interactions, each site is 
singly occupied in the ground state. A hole excitation can 
propagate freely so that the dispersion relation is given by 
Wh(fc) = —fi + 2t cos(fc). Likewise, a doubly occupied site 
can also move freely through the system. Since either of 
the two bosons of the doubly occupied site can tunnel to 
its neighbouring sites, the dispersion relation is given by 
Wp(fc) = U — fx — At cos(fc). These expressions for the quasi- 
particle dispersions are exact to leading and first order in 
strong-coupling perturbation theory. 

In fig. [7] we show the spectral functions in the 'super- 
fluid' phase for p = 1 close to the Mott transition, 
t = 0.305. The elementary excitations concentrate around 
(k = 0,oj = 0). This confirms the formation of a 'con- 
densate', as also seen in the momentum distribution, see 
fig. Moreover, it shows that the low-energy excitations 
near k = indeed dominate the spectral functions. We 
used this concept for the analysis of the ground-state cor- 
relation functions, see eq. (0]). Note that deep inside of 
the Mott phase the system size dependence of the spec- 
tral functions is insignificant, see inset of fig. El 

Conclusions. — In this work we have investigated 
the one-dimensional constrained Bose-Hubbard model 
{rib < 5) at zero temperature. Using the density- 
matrix renormalisation group method we have obtained 
the Tomonaga-Luttinger parameter from the density- 
density correlation function and determined the critical 
couplings (t/U) c = 0.305(1) for density p = N/L = 1 and 
(t/U) c = 0.180(1) for density p — 2 which separate the 
'superfluid' and Mott insulating phases. 

In the 'superfluid' phase, the momentum distribution 
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Fig. 7: Intensity (left panels) and line-shape (right panels) of 
the single-boson spectral functions A(q,u}) for t c = 0.305 ('su- 
perfluid' phase) with system size L = 64 at filling p = 1 using 
the DDMRG technique with open boundary conditions for a 
broadening r\ — 0.04 



diverges for small momenta, n(\k\ — > 0) ~ ~ L v 

(y = 1 — Kb/2), and the spectral function is finite only 
for small frequencies and momenta. In the presence of a 
confining potential, we recover the Mott plateau in the 
particle density for filling p = 0.40 and the wedding-cake 
structure for filling p = 0.54. 

In the Mott insulator, the momentum distribution is 
a continuous function. The spectral function is well de- 
scribed in terms of free quasi-hole and quasi-particle ex- 
citations which have a very long life-time for strong cor- 
relations. Their dispersion relation can be obtained from 
strong-coupling perturbation theory. A calculation of the 
quasi-particle bands beyond first order remains to be done. 



The authors would like to thank H. Frahm and T. Gia- 
marchi for valuable discussions. SE and HF acknowledge 
funding by the DFG through grant SFB 652. 

REFERENCES 

[1] Greiner M., Mandel O., Esslinger T., Hansch T. 

W. and Bloch L, Nature, 415 (2002) 39. 
[2] Bloch I., Dalibard J. and Zwerger W., Rev. Mod. 

Phys., 80 (2008) 885. 
[3] Fisher M. P. A., et al, Phys. Rev. B, 40 (1989) 546. 
[4] Teichmann N., et at, Phys. Rev. B, 79 (2009) 224515. 
[5] Freericks J. K., et al, Phys. Rev. A, 79 (2009) 053631. 
[6] Capogrosso-Sansone B., Prokof'ev N. V. and Svis- 

tunov, B. V., Phys. Rev. B, 75 (2007) 134302. 
[7] Kato Y., Zhou Q., Kawashima N. and Trivedi N., 

Nature Physics, 4 (2008) 617. 
[8] Capogrosso-Sansone B., et al, Phys. Rev. A, 77 

(2008) 015602. 

[9] Stoferle T., et al, Phys. Rev. Lett, 92 (2004) 130403. 
[10] MORA C. and Castin Y., Phys. Rev. A, 67 (2003) 
053615. 

[11] Giamarchi T., Quantum Physics in One Dimension (Ox- 
ford University Press, Oxford) 2004. 



[12 

[13 

[14 

[15 
[16 

[17 

[18 

[19 
[20 
[21 
[22 
[23 
[24' 

[25 

[26; 
[' 2 ~ 

[28 
[29 

[30 
[31 

[32 

[33 

[34; 

[35 

[36 
[37 

[38; 

[39 

[40; 

[41 
[42 



White S. R., Phys. Rev. Lett, 69 (1992) 2863; Phys. 
Rev. B, 48 (1993) 10345. 

Kuhner T. D. and Monien H., Phys. Rev. B, 58 (1998) 
R14741. 

Kuhner T. D., White S. R. and Monien H., Phys. Rev. 
B, 61 (2000) 12474. 

Jeckelmann E., Phys. Rev. B, 66 (2002) 045114. 

Lieb E. H. and Liniger W., Phys. Rev., 130 (1963) 1605; 

Lieb E. H., Phys. Rev., 130 (1963) 1616. 

Lieb E., Schultz T. and Mattis D., Ann. Phys. (N. Y.), 

16 (1961) 407. 

Frahm H. and Schadschneider A., The Hubbard Model, 
edited by D. Baeriswyl, et al. (Plenum Press, New 
York) 1995, p. 21. 

Frahm H. and Schadschneider A., J. Phys. A: Math. 
Gen., 26 (1993) 1463. 

Freericks J. K. and Monien, H., Phys. Rev. B, 53 
(1996) 2691. 

Elster N. and Monien H., Phys. Rev. B, 59 (1999) 
12184. 

Koller W. and N. Dupuis, J. Phys.: Condens. Matter, 
18 (2006) 9525. 

Batrouni G. G. and Scalettar R. T., Phys. Rev. B, 
46 (1992) 9051. 

Kashurnikov V. A., Krasavin A. V. and Svistunov 
B. V., Pis'ma Zh. Eksp. Tear. Fiz., 64 (1996) 92 [JETP 
Lett, 64 (1996) 99]. 

Sengupta P., Sandvik A. W. and Campbell D. K., 
Phys. Rev. B, 65 (2002) 155113. 

Ejima S., Gebhard F. and Nishimoto S., Europhys. 
Lett, 70 (2005) 492. 

Haldane F. D. M., Phys. Rev. Lett, 47 (1981) 1840. 

Giamarchi T., Phys. Rev. B, 46 (1992) 342. 

Moseley C, Fialko O. and Ziegler K., Ann. Phys. 

(Berlin), 17 (2008) 561. 

Shashi A., et al, arXiv: 1010.2268. 

Zakrzewski J. and Delande D., AIP Conference Pro- 
ceedings, 1076 (2008) 292. 

Ejima S. and Nishimoto S., Phys. Rev. Lett, 99 (2007) 
216403. 

Lauchli A. M. andKoLLATH C, J. Stat. Mech., P05018 
(2008) . 

Kashurnikov V. A. and Svistunov B. V., Phys. Rev. 
B, 53 (1996) 11776. 

Rombouts S. M. A., Van Houcke K. and Pollet L., 

Phys. Rev. Lett, 96 (2006) 180603. 

Kollath C, et al., Phys. Rev. A, 69 (2004) 031601(R). 

Damski B. and Zakrzewski J., Phys. Rev. A, 74 (2006) 

043609. 

Batrouni G. G., et al., Phys. Rev. Lett, 89 (2002) 
117203. 

Menotti C. and Trivedi N., Phys. Rev. B, 77 (2008) 
235120. 

Van Oosten D., van der Straten P. and Stoof H. 
T. C, Phys. Rev. A, 63 (2001) 053601. 
Pippan P., Evertz H. G. and Hohenadler M., Phys. 
Rev. A, 80 (2009) 033612. 

Jeckelmann E. and Fehske H., Riv. Nuovo Cimento, 
30 (2007) 259. 



p-6 



